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ABSTRACT 


The failure point of a unidirectional composite subjected to longi- 
tudinal shear loading is calculated. The method of analysis is based on 
the finite element technique, the incremental plasticity relations of 
Prandtl-Reuss, and the von Mises yield criterion. The failure point for 
a boron-aluminum composite and a boron-epoxy composite were comeniced to 


determine the effect of the matrix material on the composite. 
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I INTRODUCTION 


The need for material possessing a higher strength to weight ratio 
has become abundantly clear to modern industry. Existing materials 
satisfy the present high strength-low weight requirements (high per- 
formance), but composite material technology, still in its infancy 
promises to raise the standards of high performance even higher. Com- 
posite material technology is particularly interesting to the aerospace 
industry which is constantly striving to improve the performance of aero- 
nautical and aerospace vehicles. Therefore a possible solution to the 
aerospace industries’ quest for higher performance materials is afforded 
by composite materials. A composite material combines the high strength 
and high modulus of one material and the low-strength and low modulus 
of another material (or materials) to yield a material exhibiting the 
necessary high strength-low weight ratio desaned by the aerospace 
industry. 

A particularly useful composite is the two material fiber reinforced 
composite which consists of a filament embedded in a compatible matrix 
material. The filament possesses the high strength and high modulus 
properties required while the matrix, having a low strength and low mod- 
ulus adds the ductility necessary in a high strength low weight ratio 
composite. The matrix also provides the form of the composite required 
for a particular structural use. A composite plate structure usually 
consists of many layers of composite lamina with the adjacent lamina 


having different filament orientations. 
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The composite structure is analyzed considering each lamina as an 
anisotropic homogeneous material, which results in the calculation of 
macrostresses. If a lamina is examined on an individual point by point 
basis, and considering the difference in the mechanical properties of 
the filament and the matrix, the stresses calculated are microscopic. 
The microscopic stresses are important in determining at which point in 
the composite that yield or failure occur. 

To determine the yield or failure point of the composite, tensile 
tests are performed on the matrix and filament materials until the 
elastic and plastic limits are reached. The value of the limits, 
namely the yield stress and failure stress are microscopic stresses. 

To accurately determine the behavior of a composite material the macro- 
scopic stress state must be related to the microscopic stress state. 
This may be accomplished using various analytical approaches. In this 
study the finite element technique is employed. 

Several authors have considered the macrostress analysis of uni- 
directionally reinforced composites. Bloom and Wilson [4] have calculated 
the elastic microstresses for longitudinal loading. Longitudinal shear 
loading analysis was obtained by Adams and Doner [1], Baker and Foye [3] 
and Tsai, Adams, and Doner [15]. The transverse normal loading problem 
was considered by Foye [8], as well as Adams and Doner [2]. Lin, 
Salinas, and Ito [11] employed the finite element technique to solve the 
problem of combined loading, calculating an initial yield surface and 
later extending their work into the plastic region [10]. 

The present investigation considers the effect of longitudinal 
shear loading on a boron-aluminum composite and a boron-epoxy composite. 


The objective of this analysis is to determine the failure loads of the 
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two composites and to compare the effect that the two different matrix 


materials have on the composite's properties. 
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De DESCRERMONRORAREASIICIMY, 


This investigation considers the elastic-plastic analysis of a uni- 
directional reinforced composite material subjected to longitudinal 
shear loading. Elastic theory being a basic subject and plastic theory 
being a more difficult, less studied topic will be considered here and 
amplified in section three, while the elastic theory will only be 
briefly discussed in section three. This section describes some of the 
assumptions of plasticity while also describing the yield criterion 
employed in this analysis, the Prandtl-Reuss equations and the universal 


stress-strain law assumption. 


A. YIELD CRITERION 

Plasticity may be defined as that property that allows a materia! 
to deform continuously and permanently without rupture during the 
application of stresses that exceed those producing yield. The final 
distortion does not only depend on the final state of stress, but upon 
the sequence of stress states producing that final state. This is in 
contrast to elastic deformation which depends only on the final stress 
state. . 

An important assumption of the theory of plasticity of metals, 
based on experimental observations, is that the volume of material 
remains constant under plastic deformation, that is, a hydrostatic 
pressure does not cause yielding. This implies plastic strain is associ- 
ated with distortional energy only and hence the hydrostatic component 
of a state of stress does not influence the point at which yielding 


occurs. What is needed then, is some criterion, a so called yield 
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criterion, that will predict yielding under a state of stress, given 
only the yield stress determined from a simple uniaxial tension test. 

For a homogeneous, isotropic material (same properties at all points 
and all directions), a criterion should not allow hydrostatic stress to 
influence yielding. Yield must also be independent of the directions 
of the axes chosen to define the system. Therefore since the material 
is isotropic, yielding will be related only to the intensity of stress 
and must be a function of the invariant of the stress tensor. The 
invariants are the coefficients of the cubic equation, which is associ- 
ated with determining the principal stresses of the stress state. If 
yielding is unaffected by hydrostatic stress (P) then it depends only 
on the deviatoric stress. This means that the criterion must be a func- 
tion of the deviatoric invariants. 

Before proceeding with the yield criterion discussion, some terms 
need explaining. Hydrostatic stress, P, is a stress that acts equally 
in all directions and will cause only a change of volume, recoverable 
on removal of the stress. This is related to the dilatation, A which 
is defined as the change of volume per unit volume. ahedeuieenie 
stress components, Oia. are obtained as the difference of the dilatation 


from the total stress tensor, that is 


Cn re were 
or 
yl) . 
where P = Be and O44 is the Kronecker delta. 


Upon expansion, equation (2-1) gives 
o..* =o .-P (2-1a) 
Oo ep help (2-1b) 
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zz 2z 
Oxy = Oxy (2-1d) 
Os = O72 (2-le) 
CGR Gay (2-1f ) 
where the hydrostatic stress, P, may be expanded to 
a“ 
P = eek Nhase Cora Too) s (2-2) 


where the Einstein summation convention is followed, where a repeated 
subscript means summation with respect to the subscript. In this study 
_ index notation as well as cartesian notation will be exchanged freely 
instead of the indices having the usual numerical values. The devi- 
atoric stress invariants are also expressed in this manner. 


J “=o 957 xx yy o =0 (2-3a) 
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Various yield criterion exist that were verified experimentally and 

that are applicable to the problem. The von Mises yield criterion is 
the simplest and most convenient to use. It may be expressed as a 
function of the deviatoric stress invariants 

F(J,*)=0 
or specifically we will use 

Jo* = k? (2-4) 
where k is a constant. When J»~<k? there is elastic behavior, that is 


there is no change in the plastic strain. When J,~=k? there can be 


elastic or plastic strain depending upon the rate of change of Jo”. 


14 














HOT RD Ge I. Tem T Om Aaggeems 
(ots) ~~ : 
(sf-S) 
(71-8) 
ot -batingiote od ysm ,9 ,z2e1se s1sateorbuh 
(S-S) eM y.04 rigal tb q 7 


- 


betssqot 6 ora ,bewol (ot af natinevios moti smmye ntydanta, | 
yout, atid ol) .dqtuatdue cds of Jovqnes dtiw nol pam ae ansem: 

ufaart-tepnadaxe ad Iliw noltston netestin 26:T tow e8 8 

-lyveb of? ..zeufev (60! vomm fouzw oft ontyat 2eotbat <A?) 


rontem ett al tezesrqee oclé-o1 zdnsisevnt ee 


(sE-5) a aaa ll 
‘ gn f.12 ahs =< 4. - 


 * 7 o « * jf + “& {—. 
Pap ste h 


(3€-$) Le eyPee?” en ‘poo » 

wes Vi Tosnamiyaqie bertiney evew San? Yatas motisstsd Bfaly ue 
ef netiadivs bielhy zoz?h nov si? \eoldotg off ap atdapliqqee 
eerbpeesiqxe od yon 1] =. Se OF SnelnevnoD Pagw Bre Tesi 
276 evar ieert2 ot obatwe of? Yo not 
oe(> ta. 


oz Thtw ow Vfl sortase 


(4-3) > ° eu 
si durtt ,volverned’otyesfs 2! ovis Sie" eh OW Serenades 6 PW ee 
oJ nas oT ‘ye -L cotW .ateove otseslq- att af sunatd of 2? sient 


AFNo ats7 ott noqu gntiaggel aise oseald mg oteeere 


| 


Other derivable forms of Jo“ follow in terms of the principal stresses 
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205 “=o, “2t09~2t+03 ~2= 2k2 (2-5a) 
or 
209 °=3[ (01-02)?+(09~03)2+ (03-01)? ]=2k? (2-5b) 
also 
Bal i 2 Ne 2 © 2 EL EO) 2729/2 e 
205 Alloa, Ovy) or T,) Ge, ca) J+2Lo,, TOV, fo, ]=2k (2-6) 


The constant k may easily be solved for in terms of Y, the yield 
stress in simple tension, by substituting into equation (2-5b), o,=Y 
and o5=03=0 

0} °+0, 2=6k2 
or 
2Y2=6k2 and a (2-7) 
ven 
Now consider the meaning of k in pure shear. Pure shear is equiv- 


alien) to ‘the’ stress state Ts 6,>0 and o3=-0, with oo=0. Substituting 


ax 
into equation (2-5b) again 


24. 24 
(crane Cra 


and (2-8) 


i Des(a2 
( is a) 6k 


Hence k is the shear stress, at which a point in pure shear, will yield. 
The von Mises criterion implies therefore, that yielding is not 
dependent upon any particular stress component, equation (2-5b) shows 
that it depends upon a function of distortion since the differences of 
principal stresses are proportional to the maximum shear stresses. 
It is convenient to state yet another variation of the yield cri- 


terion. From equations (2-5b) and (2-7), we obtain 
1 


i 


Y= {(0}-09)2#(op-03)?+(03-0}) 21 o (2-9) 
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The scaler quantity on the right hand side may be defined as an 
equivalent stress, a, so that when o=Y, yielding occurs. Thus, any 
value of c less than Y will give rise to elastic states of stress only. 
When s=Y there may be either elastic and or plastic components, depending 


upon the sign of do. 


B. FLOW RULE 

Having selected a yield criterion it is necessary to establish a 
relationship between the stress components and the corresponding deform- 
ations. This is usually referred to as the flow rule. 

St. Venant proposed that the directions of the increments of prin- 
cipal plastic strains correspond to the directions of the principal 
stresses. Later Levy proposed a relationship to connect the plastic 


strain increment to the stress in the form 





P P P P P P P P P 

deyy “yy devy ery dez dey, XY. deve “2x 
d,=———_—_ = —————— = —— =§ —— =¥ ——— = (2-10) 

xx” Cyy ay Viz OZ Ze Xx OXY ROW BY ZiolnD LeZX 


where the superscript P denotes plastic. 
The factor da establishes the relation between quantities of the 
same differential order. From the first and third equations of (2-10), 


de P ei =da(o 


aay, eaicyy Ne (2-11) 


: P Pee 
)s de de, dao, 05, 


Adding the two equations of (2-11) and recalling that the volume change 
P P P Ris 


due to plastic strain is zero, that is, dV FEE A IEG ee On, 
P Pay ccna ; 
dey tes Soe ay we obtain 
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Similar equations may be obtained for the other components by cyclic 


permutation of (2-12). 


These equations may also be written in deviatoric form, using 


CS - 4 wl 
Pes (0, .toyytoz,) “Bho 2 yyt77) 
aoe 3 

xy oxy 


and comparing with equation (2-12), we obtain 


P = a 
de, = Moy. 
P = t a 
dey = dro, 
or, in general form 
Paves 2 5 


In the above equations, a rigid plastic material is assumed because 
no mention is made of any elastic strains. For a complete statement, 


the equation should be written 
Pe y 
(laa 2 Ob\eia a (2-14) 


where the superscript P denotes the increment of plastic strain. The 


total strain increment becomes 


(2-15) 


E is the elastic strain increment. Now, according to the 


where de5 
present argument, (2-14) refers to increments of the strains as we 

follow the deformation path of the material; we must therefore add the 
corresponding elastic strain increments - that is, although the plastic 


strain increments are, by definition, proportional to the current values 
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of the deviatoric stress O45? the contribution of the elastic strain will 
be those corresponding to the change in the stress state as we go from 
two closely neighboring states of stress. 
The complete statement of the stress-strain relation becomes 

des; a re) 
Equation (2-16) is the Prandtl-Reuss equation for the deviatoric strain 
increment. The equation derives its name from Prandtl who developed the 
elastic and plastic strain G and a equations for plane stress and 
Reuss who generalized the equation, but it was Lévy who originally de- 
veloped the plastic strain equation. 


The rule of plastic flow in terms of principal stresses and strains 


also gives 


Rey cya eee nv cidn (DAT) 











C. THE UNIVERSAL STRESS-STRAIN LAW 

The Universal stress-strain law is an assumption that the plastic 
uniaxial stress-strain relation Chae f(de,P) is the ene relation be- 
tween equivalent stress and strain, that is & = f(der) for the general 
stress state. Employing the universal stress-strain law equation (2-17) 
may be expressed in terms of equivalent (effective) strain and stress as 


ae 
d\ =C-= (2-18) 


The constant C depends upon the definitions of equivalent strain and 


stress. The equivalent strain increment may be expressed 
ap ie pa2 2 2 
de 3 llder') + (der') + (des') ? 
2) P cS 
[5 te, ; de, 32 | (2-19) 
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The numerical value of [= has been chosen so that in uniaxial stress, 

Pca ND eae ii apres : 
(where devy = deo, 5 deve we have de der. The assumption of a 
universal stress-strain law then gives a value of C=3/,. The actual 


proof of the numerical value of C may be found in Appendix B. 
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III DERIVATION OF THE STRESS-STRAIN MATRIX 


The purpose of this section is to give an explicit expression for 
roPy the plastic stress-strain matrix for the von Mises material. The 
expression obtained takes a form that can be accommodated to the finite 
element analysis. ; 

The following derivations follow closely those presented by Yamada, 
Yoshimura and Sakurai [16]. The well known equations of elasticity, 
commonly called Hooke's law may be expressed in matrix form for an 


isotropic material as 
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(3-1) 
where E is Young's modulus of elasticity and v is Poisson's ratio. The 
matrix equation may be expressed as 

{o}=E[D']{e} or since E=2(1+v)/G 
where G is the shear modulus of elasticity 
{o}=2(14v) GD" Jie} (3-2) 
where al is referred to as the elastic stress-strain matrix, {o} as 


the column matrix of stress and {ce} as the column matrix of strain. 
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Hooke's law for the longitudinal shear loading case (i.e., all strains 

equal to zero except Evz and By) greatly simpifies equation (3-1) to 
“yz Ee “yz (3-3) 
or I+} 0 5 Ze oy 

It is interesting to note that in the elastic range xz and yz behavior 

are independent. This is illustrated by equation (3-3) where the 

equations for Sy7 and Ce uncouple. 


The inverse of Hooke's law, equation (3-1), may be expressed in index 


notation as 


O.: (0) 
ein ol 6. (ee) (3-4) 
Ih Ae td 3E 
where ej are the components of the strain tensor. The 045 are the 


components of deviatoric stress which are defined as 
O45 = 94575 845 °kk (3-5) 
and explicitly expressed by equation (2-1). 
Before continuing with the derivation oF the stress-strain matrix 


and deviatoric strain components, £45 must be explicitly stated, 


Exy = ec (3-6a) 
yy = Sin (3-6b) 
E55 = E557 (3-6c) 
Exy = exy (3-6d) 
Eye = evs (3-6e) 
Bag = Bax (3-6f ) 


where the dilatation strain is 
e= fe teyy S22) 3 (3-69) 
Equations (3-6) will help better understand the Prandtl-Reuss equations 


for the deviatoric strain increment expressed in index notation as 
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des5 = O45 d\+ J (B27) 
2G 
where according to equation (2-18) 
P 
dx = cE 
Go 
: 5 Oe 
Letting H* =——, dda becomes 
det ds 
dy Cau (3-8) 
oH 
and allowing the constant C to have the value three halves 
dye Ses (3-9) 
oH 


The equivalent plastic strain increment des and the equivalent stress 


G are equal to 


ae 1 
P P P\2 
de = (3 de; dei, Na (3-10) 
ele 
a= [5(oq3 0,2)] (3-11) 
do 
H* = —p shown in equation (3-8) is actually the slope of the stress 
de 


versus plastic strain curve in the non-linear plastic region (increment 
of equivalent stress (da)/increment of equivalent plastic strain ay). 
This is shown by figure 1, where for the plastic problem the transition 
must be made from the total equivalent stress-strain curve on this 
universal stress-strain curve to the equivalent stress-strain curve for 
the non-linear range. 

For the derivation of the stress-strain matrix for this elastic- 
plastic analysis of the longitudinal shear problem the equivalent stress 
becomes 

a = Boy. to2)} (3-12) 


and the plastic strain increment as becomes 


1 
Pie eu 2 2 2 
de = {3 (dep) tde Ps )} (3-13) 
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FIGURE 1. SLOPE OF THE NON-LINEAR STRESS-STRAIN CURVE. 


The von Mises yield criterion is employed in the plastic analysis 


of the problem. It is expressed here in terms of the deviatoric stress 


invariant 
Fo § a ieee 
20° = 045 945 = 3° (3-14) 
and in differential form 
Ade.2 @ Sale 3-15 
045 O55 = 30d6 (3-15) 


Solving equation (3-9) for d& and substituting into equation (3-15) we 


have 
a ai Be cal 
eee = aan dd (3-16) 
Algebraically eliminating dosj from equations (3-7) and (3-16) gives 


2 O 2 = hao 2 A 
260, 5 (dos 04544) 90-H“da (3-17) 


Using equation (3-11) and solving for da 


05; Clana 
d, = ——————_ (3-18) 
$02 (—#1) 
3G 
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For convenience let 


S = 252(14H-/3G) 


thus 
dd = O55 desi 
(3-19) 
S 
The term Gas leas may be simplified to O56 dease The proof is not very 


long and will be presented here rather than detailing it to an appendix. 
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Now the deviatoric stress increment do. § from equations (3-16), 


(3-17), and (3-19) is equal to 


de 
262 fev rat celal tare 
do. § G(de os on ) 


AE Sh aitene, Hs eae) (3-20) 


3 


which was obtained by recalling equations (3-6) which may be expressed 


in incremental index notation as 


Go aa | td) cakk (3-21) 
3 





and substituting dess and dX into equation (3-7) and solving for do. 5 
equation (3-20) was derived. 


By definition the total stress increment do. is 








(3-22) 


Finally substituting (3-20) into (3-22) 


d 
a= CAC deck mee Onl erm oneoklcs =2 
do; 5 er 8 ideoo45 a el |) (a2) 





which may be expressed in matrix notation as 

{do} = E[D?]tde} = 2(14+v)G[D Itde} (3-24) 
Mention should be made here of the particular notation for the elastic 
and plastic stress-strain equations. The elastic equation is expressed 
in the usual tensor form while the plastic equation is expressed in 


differential form. The plastic problem is point dependent, while the 
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elastic problem is not. The value of (Py is a function of what point 


on the plastic stress-strain curve the calculation is being made while 


the value of roy is a constant anywhere on the elastic curve. 


The elastic and plastic equations are similar enough to allow the 


replacement of tory for a yielded material by oak Explicitly [0°] 
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(3-25) 
For the longitudinal shear loading case coPy becomes 
2 iS 
deere yz Oy zaZx 
p v 2 S S 
= ie (3-26) 
l+v = Fy 78> 1. On 
S S 
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The plastic stress-strain matrix Por a. unlike the elastic stress-strain 
matrix Dee does not uncouple to give independent xz and yz behavior. 
The matrix roPy is the non-linear, point dependent, plastic stress- 
strain matrix relating incremental stress and strain. 

The next section describes how [oF] and na fit into the scheme of 


the analysis to determine the failure point of the composite. 


27 





-2esue, sleet. peer intog ae ets bade 

. tarda, bon agonde Fsaeonoro an ; 
10 snort, oad. otnt 2? [07 ons (30) wort: zediaaz8b ¢ 
sortzoqmss att to tatog owiist odd pee ray 


+H - ri 


NS 






TV METHOD OF ANALYSIS 


The objective of this analysis is to determine the failure point of 
a unidirectional fiber reinforced composite material subjected to longi- 
tudinal shear loading and to determine the effect that the matrix 


material has on the overall properties of the composite. 


A. THE NON-LINEAR PLASTIC PROBLEM 
The method for the plastic analysis described here was developed 
by Yamada, Yoshimura and Sakurai [16] for plane stress and modified for 
the longitudinal shear loading case. Briefly the analysis solves the 
incremental force-displacement equilibrium equation 
[KP] {dw} = {dF} (4-1) 
where cK] is the stiffness matrix, {dw} is the differential displacement 
in the Z-direction and {dF} is the differential force vector in the Z- 
direction. Equation (4-1) is obtained by finite element discretization 
of the displacement vector {w}, this may be written 
KP (w)](dw} = {dF} (4-2) 
and solving for {dw} 
{dw} = [KP (w)]” {dF} (4-3) 
In the above problem the differential displacement is calculated on the 
basis of the preceding displacement, that is, the non-linear problem 
is solved by incremental theory, 
{aw.} = Kernen on wey (4-4) 
The stiffness matrix cKPy in the preceding discussion is a function 
of the displacement vecter {w} in the plastic analysis only. The stiff- 


ness matrix for the elastic analysis is a function of two independent 
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elastic constants and the equation {dw} = [K] {dF} may be integrated 

to give {w} = rey {F}. Hence the linear elastic problem can be done 
in one analysis while the non-linear plastic problem requires repeated 
incremental analysis. That is, the non-linear problem is approx imated 


by a succession of linear problems. This is illustrated in figure 2. 






IN A (ELASTIC RANGE): F=kw 
IN B (PLASTIC RANGE): AF=kaw 


FIGURE 2. THE ELASTIC-PLASTIC PROBLEM 
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The present discussion is limited to general terms, but a more 
detailed presentation will be related in Section 5. 

The problem considered in this analysis is the longitudinal shear 
loading case. The edge x=a is given a uniform unit displacement in the 
positive Z-direction as shown in figure 5. The resulting macrostresses 


X XZ 
element method using a linear strain triangle (6 degrees of freedom). 


S . and Syz and microstresses o._ and Svs are calculated by the finite 


Felippa [6] presents a derivation of the linear strain triangle for in- 
plane behavior and the derivation of the linear strain triangles for 


longitudinal shear behavior is presented in Appendix A. 


B. DESCRIPTION OF THE FILAMENT COMPOSITE 

The microscopic stress state in a composite material depends upon 
the arrangement of the filament within the matrix and the applied loading. 
To simplify the present analysis, we consider the special case in wnich 
the filaments are arranged in a doubly periodic rectangular array. This. 
array of filaments in the matrix may be divided into basic unit cells 
of dimensions 2ax2b, to facilitate ease of calculation. - This arrange- 
ment is illustrated in figure 3. The arrangement illustrated in figure 4 
with square cell a=b=0.0052 in., and filament radius c=0.004 in., results 
ina filament volume of approximately 50 percent of the total volume. 
The symmetry of the basic unit cell allows the analysis to consider only 
One quadrant of the basic unit cell as shown in figure 4. The quadrant 
considered has dimensions axb with a filament of radius c. It is also 
assumed that the filament and matrix material are homogeneous and iso- 
tropic and that a perfect bond exists between the filaments and the 


matrix material, i.e., there is displacement continuity along the interface. 
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FIGURE 3. RECTANGULAR ARRAY OF FILAMENTS AND A BASIC UNIT. 





FIGURE 4. A QUADRANT OF THE BASIC UNIT. 
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The finite element model used is shown in figure 6. The model is 
divided into 92 elements and a total of 215 nodes or degrees of freedom. 
The 92 element model was chosen over a 54 element model and a 200 element 
model. Previous analysis of Lin et al [11] indicated that the 54 element 
model did not yield results of sufficient accuracy, while the 200 element 
model gave very accurate results, but were not commensurate with the 


amount of computer time required. 


C. DERIVATION OF THE EQUILIBRIUM EQUATION 
In the elastic region the resulting strains obey the linear Hookian 
stress-strain relation, 


(ol eplDn I ten? 
6x1 6x6 6x1 (ee) 
where {a} is the column vector of stresses, ro“y is the elastic stress- 
strain matrix and fea is the column vector of elastic strains. The 
explicit form of the symmetric matrix tof] is given by equation (3-1). 


The stiffness matrix [k] for the elastic elements may be derived from 


energy considerations. If 


then the principle of virtual work shows that, 


Ck] {w} = au/ (4-6) 


o{w} 
The particular expression for the stiffness matrix of the longitudinal 


shear problem for the elastic problem is obtained as follows: 


U = [|ftor"taeray =| J] [tor "tena (4-7) 
V V 


where we used the fact that stress is linearly proportional to strain, 


where the fo} is the transpose of {o}. Using equation (4-5) and (4-7), 
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we obtain 

u oral fesmnayecn ‘i (4-8) 
The Pe eeoreke relation may be expressed as 

ete [Biwi (4-9) 

where [B] is the coefficient matrix and {w} is the displacement vector 
in the Z-direction. Substituting this relation into the previous 
energy equation (4-8) 

U = eu GaNe t dA{w} (4-10) 


where t is the uniform thickness. From equation (4-10) the elastic 
stiffness matrix for an element becomes 


rk] = || ce)" co£y"Teyean (4-11) 
A 


where A is the area of the triangular element. Similarly for the plastic 
element, the differential stress-strain relation is 


flan @ (eon lest (4-12) 
6x1 6x6 ~=-6x1 


where the explicit form of [0°] is given by equation (3-25). In this 
case the element stiffness matrix is 

cP] =|f cel" ro°"ce}tan (4-13) 

The stiffness matrix ckey for the elastic elements and Ck] for the 
plastic elements are assembled to form [K], the stiffness matrix of 
the whole body. The system stiffness matrix relates the incremental 
load {AF} to the incremental nodal displacements {Aw} for the elastic- 
plastic analysis 
{AF} = [K] {Aw} i (4-14) 

Equation (4-14) is basically how the plastic, non-linear problem is 


approached. The non-linear problem is linearized by taking very small 
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finite steps (increments) along the plastic portion of the stress-strain 
curve. In this way the equivalent stress for each element will not 
deviate during succeeding steps from the actual: stress-strain curve to 
any significant amount. The following discussion supplies a little 


more detail to this sketch. 


D. SOLUTION TECHNIQUE 

Using the incremental load technique the elastic-plastic problem is 
approached as follows. First the elastic problem is done in one step. 
For a test shear displacement (unity), solve [K]{w} = {F}. This equation 
holds for the elastic range only. Calculate the elastic displacements 
at the nodes and then the elastic strains, stresses, and equivalent 


E E and scale all 


stresses o for each element. Determine the maximum o 
elastic values using a scale factor of r = Vicon where Y is the yield 
stress of the material resulting from a uniaxial tension test. This is 
done to bring the first element or elements to yield. Next calculate 

the [oP] and resulting ckPy for the yielded element or elements and re- 
place the previous elastic matrices [oFy and ed. with the new matrices 
in the system stiffness matrix [K]. Next a test increment faw!} of the ) 
shear displacement is chosen (where the superscript T inside the brackets 
Penotelthe test increment). The strain increment Ne and the stress 
increments os; are calculated at each element after solving for CAWDD. 
Now scale all elements remaining in the elastic state by the minimum of 


r+ (r244(ae 1)” (y2-52))2 
Rees pee (4-15) 

Aline > ))Z 
where 


Tl =M(Rae nny 2a eeocng URS 2 (4-16) 
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mae 10,2172 (4-17) 

R is the scale factor that brings the next elastic element to yield 

(see Appendix B for derivation). Also, o is the present equivalent 
stress of the elastic element and age denotes the increment of co induced 
by {aw'}. 

A fine point should be made clear before proceeding, ie and nel are 
not the same quantities. A significant difference does exist between 
the change in effective (or equivalent) stress (ial) and the effective 
change in stress ng). The change in effective stress is the increment 
of equivalent stress (c) caused by an incremental increase in the loading. 
The effective change in stress is that portion of the stress increment 
that actually produces a change in the plasticity of an element. 

At this point in the procedure the nodailldispilacenentes strain and 
stress increments are multiplied by the minimum value of R and added to 
the existing values. Next the equivalent stress o, given by equations 
(2-9) and (3-11) and the incremental equivalent plastic strain, Wee care 
calculated for each element. 

eur O45 hei; 
If the value of the dee is positive, the analysis returns to the be- 
ginning of the next plastic cycle by calculating [oP] and [kP] for the 
next yielded element. If ae is negative the computation is stopped. 
This is included here because large plastic strain are possible without 
any practical increase in the load after the plastic region has expanded 


tremendously. This gives an indication that the analysis is very close to 


the collapse load of limit analysis. The problem may terminate in a very 
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natural mode when all the cycles are complete or when the incremental 
equivalent plastic strain, na goes negative. These two modes of failure 
are illustrated in figure 7. The next section describes the actual 


procedures in much more detail. 
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FIGURE 7. MODES OF FAILURE 
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V PROCEDURES OF CALCULATION 


This section provides a much more detailed explanation of the pro- 
cedures of calculation than does the Method of Analysis section. The 
procedures of calculation of this analysis was suggested by Yamada et al 
[16] and modified for the longitudinal shear loading problem. The analysis 
was constructed to be able to perform calculations for both the xz and 
yz shear loading problems, but will only be described for the xz shear 
problem. The xz shear problem was the one actually considered by this 
investigation. The procedure follows. 

For the elastic analysis (Kw=F), 

1. Apply a uniform unit displacement to the x=a face, and calculate 
the resulting displacements, strain, and stresses for tne elastic be- 
havior (system stiffness uses elastic stress-strain relations here). 
Figure 5 illustrates the uniform unit displacement given to the basic 
unit cell (first quadrant only). 

2. Scale the elastic values so that the element with the maximum 
equivalent stress Gre) is at initial yield (Y), i.e., adjustment of 
input displacement so that a single point is at yield. The scale factor 
re is 


ES 
max 


and the nodal displacement vector {we} and the nodal load vector ip at 
the initial yield point are given by 
Ey ety) | (5-2) 

i (5-3) 


a 
= 
w 
iT] 


rE tr 


co 

= 

— 
YW 
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This terminates the elastic analysis. For the plastic analysis that 
follows we proceed in cycles or increments. The procedure for the plastic 
analysis is one cycle that is repeated over and over until one point or 
element of the composite has failed (termination of the problem may end 
in other ways, see section four and actual computer program). The loading 
for the plastic problem may mean either displacement or force loading, 
but for this particular study the displacement loading was employed. 
Displacement loading is used because we are interested in achieving the 
shear derouracion <tate. Thus we use equation (4-4) for this problem 
({dw,} = [K(w,_,)]-*{dF}). 

3. Calculate [oP] and ck] for the post yield elements; usually 
only one element will reach yield in the first cycle. However since 
elements close to yield will be carried into yield during the next incre- 
ment of displacement, any element that is within .995 of yield is con- 
sidered to be at yield. The value .995 is completely arbitrary; it 
depends on the accuracy that the user wishes to obtain. 

4, Modify the system stiffness matrix to reflect the stiffness of 
the elements that have yielded. 

5. Choose a test increment aw! of displacement. 

| 6. Solve {dF} = [K(w._,)] (dw) and then calculate the strain 
increment {Ac} and the stress increment {Ao}. 

7. Calculate the scale factor R (equation 4-15) for each remaining 
elastic element to reach yield. 


8. Determine the minimum R. The displacement increment, R faw!} 


min 
is sufficient to yield that element having the minimum R and place it in 


the plastic range. 
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9. Multiply the strain increment and stress increment calculated 
in step six by Ra and add the resulting del (A) quantities at the 
nodes to the previous nodal values, i.e., {w}i={w}, ,+{Aw}. Store the 
results. 

10. Calculate and store the new equivalent stress for each element. 

11. Calculate the increment of plastic strain for each post yield 
element. The stiffness of each plastic element constantly changes as the 


element progresses into the plastic range. 


, 


CoS 55 mens 
PPLE palate |) (4-18) 
G(1+H"/ 5) 
Bor achiiSmanallySiSie- sen (oh Net oI Ness }) 
at 2 yZ YZ ZX ZX (5-4) 
o(1+H"/,-) 


12. Check that Aue is positive. If ae. As positive return to step 
three. If Ie ae negative stop computation. This is included because 
large plastic strains are possible without any practical increase in the 
load after the plastic region has expanded tremendously. This is one 
way of ending the problem, the other way of course is a natural ending. 
A natural ending occurs when the number of plastic cycles requested are 
run or where some of the other checks built into the program terminate 


the problem (see computer program in Appendix C). 
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VI RESULTS, CONCLUSIONS AND RECOMMENDATIONS 


Aer RESUIETS 


The computer program for this elastic-plastic analysis of the 


longitudinal reinforced composite loaded in the xz shear mode was run 


until some point in the matrix reached its ultimate strength. 


The 


properties for the matrix material and for the boron filament used in 


the analysis are given in table 1. 


In the analysis the first point (element) to yield was calculated 


before the ultimate breaking point was located. 


The initial yield area 


for boron-epoxy and boron-aluminum composites occurred at element 50 


as shown in figures 8 and 9 respectively. 


Bordering the filament-matrix 


interface, this point (element 50) lies in the matrix material on the 


lower horizontal face of the model. The filament point in both com- 


posites that obtains the highest equivalent stress at the time of initial 


yield is element 32, which is the filament element adjacent to matrix 


yield element 50, (see figure 6). 


Table 1. 
Property . Boron 
Young's Modulus 60x10° psi 


Poisson's Ratio 0.20 
Elastic limit 480,000 psi 
Ultimate Strength 480,000 psi 
N = 
B a 


* Tensile Properties for Epoxy 


4) 


Material Properties 


Aluminum 
10x10° psi 
0.30 
34,600 psi 
42,500 psi 
865) 
46,300 psi 


Epoxy* 
.5x10° psi 
0.31 

3000 psi 
4200 psi 
Sy, 
5088 psi 
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When the analysis finally reaches the breaking point for both com- 
posites, some interesting results have occurred. In the boron-epoxy 
composite element 50 becomes the first matrix element to fail and 
element 32 is still the highest stressed filament element, but in the 
boron-aluminum composite the initial yield element 50 is overtaken by 
element 70 to become the first failure element and element 32 is over- 
taken by element 31 to become the highest stressed filament element at 
the time of failure. The progression from initial yield to failure for 
the boron-epoxy and boron-aluminum composites is shown in figures 8 and 9 
respectively. 

The results of the initial yield analysis are presented in table 2. 
These are the results associated with the elastic limit analysis. The 
boron-aluminum composite has a macro shear stress S\ 7713372 psi at its 
elastic limit; and the boron-epoxy composite has a macro stress of 


S, 771044 psi at its elastic limit. 


Table 2. Results at Initial Yield 


Boron-aluminum Boron-epoxy 
x7 65 SS 9535 psi 1729 psi 
Cmaxsg 34,600 psi 3000 psi 
sian 35,288 psi 3070 psi 
Sie 1S Ovica DSi 1044 psi 
Syz 240 psi 22 psi 


We note that the maximum stresses in the filament are but a small fraction 
of the elastic limit stress (480,000 psi) of the filament material, that 


is, the filament is not being used efficiently. Thus, designing to the 
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PROPAGATION OF PLASTIC ENCLAVE FOR BORON-EPOXY 
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Figure 9. PROPAGATION OF PLASTIC ENCLAVE FOR BORON-ALUMINUM 
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elastic limit of each, the boron-aluminum composite can sustain a shear 
load 1300 times greater than the boron-epoxy composite. 

Table 3 contains the post yield results. Continuing the loading 
into the plastic region, the more ductile aluminum matrix allows the 
boron-aluminum composites' longitudinal shear loading, Sie to increase 
80 percent (from 13,372 psi to 24,140 psi), while the boron-epoxy 
composites' Syy increases 50 percent above its elastic limit (from 1044 psi 
_to 1539 psi). These increases in longitudinal shear load correspond to 
a 25 percent increase in the equivalent stress of boron-aluminum (the 
equivalent stress in the aluminum element 50 goes from 34,600 psi to 
42,500 psi), while boron-epoxy's equivalent stress increases 40 percent 
(o increases to 4200 psi from 3000 psi). Clearly the boron-aluminum 
offers the designer a better shear load carrying ability than does the 


boron-epoxy composite when the design is allowed to reach failure. 


Table 3. Results at Ultimate Strength 


Boron-aluminum Boron-epoxy 
x76 5 24,444 psi 2427 psi 
x70 24,563 psi - 
Gnaxso 42,372 psi 4212 psi 
Cmax 42,547 psi - 
Cmax 51,097 psi 4356 psi 
Onaxa4 52,732 psi - 

Ses 24,140 psi 1539 psi 
Syz 692 psi 33 psi 
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An all aluminum unit cell was subjected to the longitudinal shear 
loading to verify the analysis' reliability. A Sy macro shear of 
20,000 psi was obtained for an equivalent micro stress, co of 34,600 psi 
at initial yield. The homogeneous cell had a uniform microstress state 
throughout, with o,. = o/J3 = 20,000 psi and Sip =e The ultimate 
macro shear stress Sso54 is 24,500 psi with a uniform equivalent micro 
stress o = 42,500 psi, resulting from ie 24,000 psi and Ove = 0 psi. 
We note in the case of a homogeneous cell subject to pure xz shear, 
there results only xz behavior, that is, in this special case post yield 
behavior does not give coupling between xz and yz behavior. This is 
evidenced by Ove = 0 psi. 

The composite results were further verified by constructing an 
equivalent stress-strain curve (universal stress-strain curve) and then 
by plotting the values of equivalent stress and strain from the computer 
program on the same figure to show how close the analysis came to the 
input stress-strain curve. The elastic portion of the curve was con- 
structed using the linear elastic relation, G = Ee, while the plastic 
portion of the curve was constructed using the power law, (oye 
The portion of the boron-aluminum equivalent stress-strain curve where 
the slope becomes constant (i.e., for & > 40,500 psi) was plotted using 
the relation, H- = cor The power law values, N and B were calculated 
using the data pended original uniaxial stress-strain curve for the 
material. The determination of N and B for the aluminum was easily ob- 
tained, but N and B for the epoxy curve became extremely difficult to 
calculate due to the almost linear stress-strain curve for epoxy in the 


plastic region. Figures 10 and 11 show the equivalent stress versus 


total strain curves for the boron-epoxy and the boron-aluminum 


x 
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respectively, while figures 12 and 13 show the equivalent stress versus 
the equivalent plastic strain. The close proximity of the plotted 
values to the curves verify small enough incremental steps were taken 


during the non-linear plastic analysis. 


B. CONCLUSIONS 

In this analysis the objective has been to develop and use a program 
which could determine the failure point of a unidirectional reinforced 
composite material subjected to jong atudinad shear loading using the 
finite element method and to make a comparison of two composites. 

The failure point was determined for the two composites and occurred 
at element 50 for the boron-epoxy composite and at element 70 for the 
boron-aluminum composite. A greater portion of the aluminum matrix 
yielded before failure than for the epoxy matrix. This occurrence came 
as a result of the more ductile aluminum matrix. 

The appearance of the longitudinal macro shear Syz is smal] compared 
to the longitudinal shear Syo? but must be accounted for because the 
Oye microstresses were comparable to the oe microstresses for the inter- 
face elements (values of 50 percent or greater were reached). Since the 
significant Svz stresses occur in the interface elements when the model 
was loaded in Se shear ,xz- yZ coupling may be explained as a result of 
the circular filament cross section. If the filament cross section had 
been rectangular the Ov7 terms and hence Syz term would have been less 
likely to appear. 

Considering all the results the boron-aluminum composite definitely 
seems to possess a better longitudinal shear load carrying ability than 


does the boron-epoxy composite. But the all aluminum cell gives the 


same results. Hence it would not be economically wise to use a boron- 
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aluminum laminar composite in place of a homogeneous aluminum plate if 
in-plane shear loading was the only load applied. The boron-aluminum 
composite may not fare to well when compared to an all aluminum material 
in the longitudinal shear mode, but from the previous analysis of Lin 

et al [10] we see that the all aluminum material would not approach the 
boron-aluminum composite when loaded in other modes (i.e., transverse and 


longitudinal loadings). 


C. RECOMMENDATIONS 

The logical extension of this analysis would be to combine it with 
the other loading problems studied by Lin et al [10] to construct a fail- 
ure surface. Also the investigation of different shaped cross sections 
for the filament should be considered. These suggestions combined with 


a parametric analysis could lead to other interesting thesis topics. 
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APPENDIX A 


CALCULATION OF THE STIFFNESS MATRIX 
USING FINITE ELEMENTS 


Section three derived the stress-strain matrix for the plastic analysis 
and stated the stress-strain matrix for the elastic analysis. This sec- 
tion will show how these stress-strain matrices are expressed in finite 
element form so that they may be written in computer language. 

In this analysis, the longitudinal shear loading problem is associated 
with the macroscopic stresses Spas and Syz and the displacement u=v=o, 
w=w(x,y). The finite element used in this analysis is the linear strain 
triangle (LST). For the longitudinal shear element the elastic strain 


energy is given by 


U = 1) 0452450" (A-1) 
V 


which becomes for this problem 


5 il 
u 5 || (22,28 242%yzt yz AV 


. V 
meee 
v = EI) oz txz*yatyz)0¥ (A-2) 
V 
where Gos and Ove are the microstresses and lie and Yyz are the shear 


strains associated with the problem. Assuming a unit thickness 
auk. ey = 
5 |] 22 ya) xz dA 
7 
yz 


A 
yy x2 [T\Yx2{daet | | fo}! ty 3A (A-3 
2 Oo y 

A yz yz A 
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The minimization of the strain energy U with respect to the displace- 
ment vecter {w} gives the stiffness matrix of an element. The derivation 
of the element stiffness matrix follows. 


In the elastic elements, 


Ore Gy, and Oy2 = Gry, 
oO Gn iy, 
a 4 pte een ae 
yz yz 


Substituting into equation (A-3) 


u - feat (A-5) 


Before the derivation may continue the array {y} must be defined. 


ga Teale | Se ee YXZy 
: 
y2 1s Ba TXZ5 
; Kee 
XZ, 
Vyzy 
‘yz, 
i 
YZ 
y 
= (a ‘ XZ ; 
Dee an, 23 3 (A-6) 
2x 6 ‘yz, a 


where vy, are triangular coordinates. 
The Neos and Yyz $ go from one to three for the corner nodal points 
of a linear strain triangle. Now the nodal shear strains may be ex- 


pressed in terms of the nodal displacements 
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where the Ny and Ny are the first partial of the shape functions with 
j i 
respect to x and y respectively. Evaluating equations (A-7) and (A-8) 


at each of the six nodes, 


WxZhy eG = Ged ko? and (A-9) 


we. bcc ft (A-10) 


Combining results and putting equations (A-9) and (A-10) in matrix form 


occas = [zx] {w} 


A-11 
6X1 6X6 6X1 ( ) 
LY yz} = [zy] {w} (A-12) 
6X1 6x6 6x1 


Substituting equations (A-11) and A-12) into (A-6) we have 


[zx] 
{y} = [é] {w} (A-13) 

[zy] 

= [e] fz] tw 
Now using the equation (A-13) and equation (A-5) we finally have 
U = £) hore ter WF" EsICzd nen (A-14) 


Letting [B] = [e][z] we have 
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u = Srew)! Cay" DoF y! 
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[B]{widA (A-15) 


Minimizing U with respect to {w} and integrating 
rk] = [B]'[o'ILB]ta (A-16) 
sffeichh is identical to equation (4-11) 
For the plastic elements the stiffness matrix becomes more involved. 
The stress-strain relatives are a function of the loading (see equation 
(3-26) ). The only difference is that [oP] takes the place of [oF] in the 


plastic stiffness matrix, thus 


ckP] = ce]'coPipejea. (A-17) 
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APPENDIX B 


SOME DERIVATIONS OF EQUATIONS 


In this appendix will be presented some important derivations necessary 
to the understanding of the elastic-plastic analysis. There are two 
proofs of the scale factor R used in section four and five. The first 
proof is geometrical in nature, while the second is an algebraic deriv- 
ation. The derivation of the gamma, which is necessary for the scale 
factor R derivation is also included. Finally it is shown that the con- 


stant C of equation (2-18) equals 3/2. 


A. GEOMETRIC DERIVATION OF EQUATION (4-15) - 

Here we obtain equation (4-15) through a geometric viewpoint. The 
basis of the discussion is contained in figure 14. In figure 14 c is the 
present equivalent stress of the elastic element and Neu denotes the in- 
crement of c induced by the displacement increment faw!}. o and ee are 
represented respectively by OP and PS, while Ao i PR. For the longi- 
tudinal shear problem (3+i0")” (see equation (B-7)) and (Ao Ty2 (see 
equation (B-13)) determines the load increment sufficient to just cause 
yield in each elastic element. 

PR= We | Re 
Using the cosine law, 


Let x = o2+(As 
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STRESS PATH OF THE ELEMENT 


FIGURE 14. DETERMINATION OF THE LOAD INCREMENT TO 
JUST YIELD THE ELASTIC ELEMENT. 


and 


COSa = x = ARCCOSx SINa = ¥ 1-x? 
Using the sine law, 


et came = 180 = (ato) 


SINo SIN@ SINO 


Qd = Y ORR= aa 


Son eo eee Pe 
/1-x2 — SIN@ 1-x2  SINo 
Sills 2 @ I a 
Y 
@ = ARCSINZ 
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¢ = 180 - (ARCCOSX + ARCSINZ) 


so Mo pO os AS GLI el ae 
J]-x2_ SIN(180-(ARCCOSX + ARCSINZ)) | SIN(ARCCOSX + ARCSINZ) 


LET g = ARCCOSX and Y = ARCSINZ 
aERO 
SINBCOSy + COSgBSINy 
_ zeal 
Ly TET) PCNA. eae 
XZ 4 
1-x2 


3 vl 1-(8) (1x) + x) fese 
Ji=x2 
[D Tey) ‘ 
Y -& (1-x ) + x6 


i 
PQ = {£(v2-62)4(as ")°4x2]/4(aa ")°3? + 





ee} 
Oo 
iH 
< 
oa 
— 
1 
N 
NS 
+ 


(B-1) 


2 2 
Gee “sexs yan } (8-2) 
Expanding terms and letting 
T = (Ac 1)? ogas! - Gisw ; (B-3) 
PQ becomes 
ay ihe T 
PQ = [{(Y2-e2)4(Ho | )2+r2}°+r]/2i6 (B-4) 
i 
8 [MWese2)a(ws ane) ae) 
ie MP (B-5) 


B. ALGEBRAIC DERIVATION OF EQUATION (4-15) 
The next derivation of R may be expressed in general terms but is 
instead derived with the particular problem of longitudinal shear loading 


in mind. 
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displacement RVine 


Fd 2 2) S 
o 3[o,, toy, ] (B-6) 


BU arco ilk 
(GtAc )2 = aie, lie, POG )2] (B-7) 
ey Sate wera oeeul| if 

2 QS 2 
62+26 Ao t+(Ac ) 3[o,, +20, Ao, 5 


2 
ate. Ao ne ey 


Ty2 
DANS Sou yz yz 


yz 
Also when ioc 2 = Y2 element is at yield: 


T mer Slt 
De 2 = 2 2 
Y 3[o,, +2Ro, AG), +R (a5, ) 


hte Stil 
2 2 2 = 
+ ove ANS) As, +R (Ao, )2] (B-8) 


e tigi fs 
2 = SP 
Y o +3[2R(o, 05,5 toy, AG, 


DEN = oatee es 2 
+R (ast thoy} )] 

Mea 
ox oxze yz acyz 


(Bo")? = 30 (a5,,")2#(05,.")?] (B-10) 


r= -6[ (B-9) 


‘Making the necessary substitutions, 


(a! )2R2-rR+82-Y2 = 0 


R = r+y r244(Ao! )2(¥2-62) i (B-11) 


2(io!)2 


C. DERIVATION OF EQUATION (4-16) FOR r 
From the original postulate 


fe (el ene a(nelye | (B-12) 
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which equals PR in figure 14. The term Xo! becomes upon expansion 


eo ae ane -q2 “72 ane alin 
WG aS elie ¥ho,,, *Ao7, *2((Ao, ) +(ho,, ) 
“T al 
(No yay )ye (B-13) 


zx 
The equivalent stress vector is 
Oma Oe ie ito. 2a to.) to. 2) 
(eens 2 = Gente del eae! 2 (B-14) 
putting o + nee in the form of the equivalent stress vector 


“72 22 
2(atas')2 = (6, 403,.") +5405") + 


XX XX shy TSAI 
(5,,#03,,') +20 (3, H06,, 1) (6, #45,,") 
i G0) 
36a8'+2(a3")? = S (aaa! )2-452 (B-15) 


Now the outcome is subtracted from (Ko")2 to give r. For this analysis 
T becomes 


ip Ss ine \e-ene! anal \2 


= -6[5. Ac. + 
6L3,. Ao, CeunNG 


y Toyz S%y2 (B-16) 


D. DERIVATION OF C IN EQUATION (2-18) 
The proof that constant "C" of equation (2-18) is simple and straight- 
forward. Equation (2-19) for a tension test and principal strains 


becomes 
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and from equation (3-11) the equivalent stress becomes 


“aponemuoeieal BOO Moen, 
o SG. O44 )2 = |3ta, too to3 2 (B-18) 


where Bi he and Ga are the principal deviatoric stresses. Equation 
(B-18) simplifies to 

= ale 

= J F ((G01)?+(-g0 )2+(-901)232 = oy 


2 


and for tension o, = <3. Now substituting into equation (2-17) 
dyeeidew: <del = Bide. 
1 Bes Ca la 
on xO Sado) 
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APPENDIX C 


THE COMPUTER PROGRAM 


The computer program performs an elastic-plastic analysis for a uni- 
directional reinforced composite material subjected to longitudinal shear 
loading. The main sections of the program are listed and described below. 

Main - Constructs the elastic stiffness matrix and modifies it according 
to loading and boundary conditions. It calls the equation solver and then 
scales the elements in the elastic and plastic analysis and controls 
further calling sequences. 

Post - Constructs the plastic stiffness matrices for the elements in 


the plastic range. 


Caltau - Calculates element strains and stresses: c=Bw. p=ee and 


Bandec - Solves banded matrix problems (Kw=F, or K(Aw)=AF). 
MULI - Multiplies a symmetric banded matrix with a vector. 
| This portion of the discussion of the computer program is devoted to 

assumptions made in writing the program and are listed here for 
convenience. 

ile The properties of the matrix material may be different in com- 
pression and in tension, YMC and YMT respectively. This assumption 
explains the dilatation variable SII's inclusion in the program. SII 
was the variable that determined if a matrix element was in compression 
or tension. It was included to preserve the generality of the program 
but does not enter the analysis in the particular case of longitudinal 


Shear loading. 
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2. If an element is within a particular fraction (called TEST) of 
the yield value, it is included as a yielded element. In this analysis, 
TEST was chosen as 0.995. Therefore for calculations any element which 
is at .995 times the yield value is at yield (TEST = .995). After an 
increment it can be determined whether an element which was close to 
yield and assumed to yield, actually did. In other words, the size of 
the incremental step and TEST should form a consistent pair. 

3. For purposes of scaling the load to take the critical element to 
yield, the scale factor Bs is restricted to be less than some maximum 
(RYMIN). This is imposed to limit the size of an incremental step taken 
along the stress-strain curve. If too large size steps are taken, the 
elements in the plastic region will provide inaccurate results (see 
figure 2). Also at a certain point in the stress-strain curve larger an 
values will be possible which would reduce the calculation time. These 
maximum values for Rae along the stress-strain curve must be altered 
according to the materials being used. They are also a result of experi- 
-ence in working with a certain material to produce better results. 

4, Also in scaling the load to take the critical element to yield 
the filament elements will be ignored. The yield of a filament is much 
higher than the yield value for a matrix element. 

5. For purposes of computer efficiency the elements of the stiffness 
matrix were derived to facilitate their use in the construction of the 
stiffness matrix in the main program and the strain column matrix in the 
subroutine Caltau. Therefore the zx and zy variable arrays listed are 
the first partial of the shape functions with respect to x and y 


respectively. 
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G2 For Hae= devdey a two point curve fitting scheme is used to 


approximate the plastic portion of the stress-strain curve. For example 
assuming the plastic strain-stress relation can be expressed as the power 


low zg? = (2) 


The values of B and N may vary for different material and for tension 
and compression of the matrix. The calculation of N and B results from 
solving the equation <= () 


curve. The selection of the two points is important and should approxi- 


for two points on the plastic portion of the 


mately divide the non-linear portion equally. The points are chosen in 
this manner so that large errors do not result. This method is an accur- 
ate and simple procedure. 

Other curve fitting expressions giving more accurate results may be 
used. It is felt that more accurate results are not necessary and do not 
yield that significant of an improvement in the results for the amount of 
time that would have to be expended. 

7. Data input is another important aspect of programming. The only 


possible trouble area exists in numbering the nodes of a triangular 
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element. The usual convention is followed by numbering the corner nodes 
first and then the intermediate nodes in a counter-clockwise fashion. 


For example, the numbering sequence should appear as: 
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